function lower_bound = kmeansBayesLowerbound(mix, data, vars)


% E[ln p(Y|C, Z, mu)]
q1 = ElnpY_CZmu(mix, data, vars);

% E[ln p(Z|pi)]
q2 = ElnpZ_pi(mix, vars);

% E[ln p(pi)]
q3 = Elnppi(mix);

% E[ln p(mu)]
q4 = Elnpmu(mix);

% E[ln q(Z)]
q5 = ElnqZ(vars);

% E[ln q(pi)]
q6 = Elnqpi(mix);

% E[ln q(mu)]
q7  = Elnqmu(mix);


lower_bound = q1 + q2 + q3 + q4 - q5 - q6 - q7;

% fprintf(' q1=%.10g \n q2=%.10g \n q3=%.10g \n q4=%.10g \n-q5=%.10g \n-q6=%.10g \n-q7=%.10g\n\n', [q1 q2 q3 q4 -q5 -q6 -q7]);
